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A scheme for separating the high- and low-frequency molecular dynamics modes in Hybrid Monte Carlo (HMC) 
simulations of gauge theories with dynamical fermions is presented. The algorithm is tested in the Schwinger 
model with Wilson fermions. 



At this conference, much discussion focussed on 
the dominant systematic error facing dynamical 
fermion simulations; the chiral extrapolation. It 
is proving to be a significant challenge to run Wil- 
son fermion simulations of QCD using the HMC 
algorithm at light quark masses 

The standard implementation of HMC in- 
troduces pseudofermion fields to mimic the 
fermion action, and generates new proposals to 
a Metropolis test by integrating the molecular- 
dynamics (MD) equations of motion for a Hamil- 
tonian in a fictitious simulation time. The maxi- 
mum step-size for a useful Metropolis acceptance 
rate is set by a characteristic time-scale of the ac- 
tion used in the Hamiltonian Recent studies 
find that as the fermion is made lighter, the force 
induced by the pseudofermions generates increas- 
ingly violent high-frequency fluctuations. This ef- 
fect is called "ultra-violet slowing-down" in Ref. 
m. In this report, we describe a modification 
to HMC that attempts to address this problem. 
A number of interesting algorithms have been de- 
vised recently that reduce the influence of the UV 
fermion modes on either the Monte Carlo update 
scheme or in the lattice discretisation 

1. MULTIPLE MD TIME-SCALES 

A scheme for integrating the equations of mo- 
tion for the MD phase of the HMC algorithm by 
introducing different time-scales for different seg- 



ments of the action was introduced in Ref. [gj. 

The leap-frog integrator is constructed from the 
two simple time-evolution operators generated by 
the the kinetic and potential energy terms. Their 
effect on the system coordinates, {p, q} are 



VT{AT):{p,q} 
Vs{AT):{p,q} 



{p, q + At p} 
{p^At dS,q} 



(1) 



The simplest reversible leap-frog integrator is 
then 



V{Ar) = Vs{^)Vt{At)Vs{^). 



(2) 



If the action (and thus the Hamiltonian) is split 
into two parts. 



n^ T{p) + Siiq) + S2{q\, 
Til 7^2 



(3) 



then the two leap-frog integrators for these two 
pieces can be constructed as 

V,{At) ^ VsA^)Vt{At)VsA^) 
V2{At) = VsAAt) 



(4) 



and a reversible integrator for the full Hamilto- 
nian can be constructed by combining the two 
schemes: 



V{At) = V2{^) 



rAr, 

Vii — ] 
m 



At. 



(5) 
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with m E Z . This compound integrator effec- 
tively introduces two evolution time-scales. At 
and Ar/m. 



We suggest that the multiple time-scale scheme 
is only helpful when two criteria are fulfilled si- 
multaneously: 

1. the force term generated by Si is cheap to 
compute compared to that of 5*2 and 

2. the split captures the high-frequency modes 
of the system in Si and the low-frequency 
modes in 5*2. 

The popular implementation of the scheme in dy- 
namical fermion simulations with pseudofermions 
is to split the Hamiltonian into the Yang-Mills 
term and the pseudofermion action: 



51 — Sg 

52 = (/ 



(6) 



Unfortunately for light fermions, the highest fre- 
quency fluctuations are in the pseudofermion ac- 
tion, which also has the more computationally ex- 
pensive force term, thus the criteria are not met. 
The key to using the method then is to implement 
a low-computational cost scheme that separates 
the high- and low-frequency fermion molecular- 
dynamics. 

2. POLYNOMIAL FILTERING 

A very low-order polynomial approximation of- 
fers a cheap means of mimicking most of the 
short-distance physics of the fermion interations. 
Recent interest in polynomial approximations to 
a matrix inverse was inspired by the multi-boson 
algorithm of Ref . [|lOj . This led to the develop- 
ment of the Polynomial HMC (PHMC) algorithm 
pT[ | . Following this idea, we write an exact repre- 
sentation of the two-flavour probability measure 



det M = / V(t)V4)* VxVx* e 



with 



(7) 



S^ - \[MV{M)]-^^\^ and S^ = \V{M)x\^ (8) 

The fields are modified pseudofermions and 
we term the new fields x, "guide" bosons. The 
algorithm exactly recovers the probability mea- 
sure for the two-flavour theory for any choice of 
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Figure 1. The HMC acceptance rate as a func- 
tion of MD step-size for the standard leap-frog 
method, and for the new algorithm with iV = 4 
and 16 



polynomial. This scheme is similar to the split- 
pscudofcrmion method of Ref. |^ with one dis- 
tinction; using polynomials allows us to modify 
the fermion modes easily and cheaply. The suc- 
cess of the algorithm hinges on the empirical ob- 
servation that the fermion modes that induce the 
high-frequency fluctuations are localised and can 
be handled by very low-order polynomials. This 
means the force term for the guide bosons, gener- 
ated by S^ can be computed cheaply and a multi- 
ple time-scale integrator can be built that simul- 
taneously fulfills the two criteria discussed earlier. 
The time-scale split is then 



5 = S, 



G 



Scb 



(9) 



UV; as cheap ir.., os costly. 



Note that in general, p time-scales can be in- 
troduced straightforwardly by first adding p — 1 
sets of guide fields, then constructing a heirarchy 
of leapfrog integrators, with rrij leap-frog steps 



at the 
At > - 



•th 

At 



level. 



This introduces p time-scales, 

Ar 



3. TESTING THE ALGORITHM 

The performance of the method is under in- 
vestigation in simulations of the two-flavour mas- 
sive Schwinger model (2d QED). HMC runs pre- 
sented here are performed on 64 x 64 lattices, with 
/3 = 4.0, K = 0.2618. Non-hcrmitian Chebyshev 
polynomial approximations are used Mea- 
surements of the dependence of the acceptance 
probability on the polynomial degree and tuning 
the MD parameters, Ar and m are being made. 

Fig. |l| shows the acceptance rate for the 
standard HMC algorithm (with even-odd pre- 
conditioned pseudofermions) along with the 
polynomial-filtered method. Two different poly- 
nomials were used; TV = 4 and TV = 16. The solid 
lines are fits to the expected behaviour of HMC 
as At 0, namely 

(Pace) =erfc{(Ar/To)'}, (10) 

where the value of tq is determined by the best fit. 
To is then a characteristic time-scale for the modes 
encapsulated in the pseudofermionic part of the 
action, in Eqn. ^ Since evaluation of the force 
term dStj, dominates in this parameter range, tq 
is a good indicator of algorithm performance. 

4. TUNING THE ALGORITHM 

The algorithm has a number of free parame- 
ters, allowing a good deal of scope to optimise 
the performance. Fig. |^ shows the characteris- 
tic time-scale for pseudofermion integration, tq 
as a function of Chebyshev polynomial degree. 
To rises very rapidly initially as the degree of the 
filter algorithm is increased, and tq for = 16 
is about four times larger than for the standard 
HMC algorithm. Inverting the pseudofermion 
matrix requires roughly the same computational 
cost, and thus the algorithm is four times more ef- 
ficient than HMC (assuming autocorrelations are 
the same). 

For low values of N (the polynomial degree) the 
computational bottleneck is solver performance 
while for larger values, the evaluation of the force 
term arising from interactions between the gauge 
fields and the guide bosons begins to dominate. 
Effective and simple strategies for tuning the al- 
gorithm are still under investigation. We are also 
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Figure 2. The characteristic time-scale tq vs poly- 
nomial degree, N. The dashed line guides the eye. 



investigating alternative choices of polynomial fil- 
ters beyond Chebyshev approximation. 
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